
rm(list=ls())

library(arm)
library(dplyr)
library(stargazer)
library(xtable)
library(margins)
library(summarytools)

######################################## LOAD DATA and SET VARIABLES ########################################

load(file="/Users/noahbuckley/Dropbox/MayorArrest/Replication data files/Final files/replication data.Rdata")
dat <- within(dat, {
  pop_margin <- logpop*mayormargin
})

# Set up parameters for latex output
omstat <- c("aic", "bic", "ll", "f", "ser")
tplc <- "htb!"
alp <- 0.7
csw <- "-20pt"
note <- c("$^{*}$ indicates p$<$0.1; $^{**}$ p$<$0.05; $^{***}$ p$<$0.01", "Note: Coefficients shown are average marginal effects from logistic regression with HC0 standard errors", "clustered at the region level.")
dv <- "DV: Mayor Arrest"
dig <- 2

######################################## DESCRIPTIVE TABLES #################################################

############# TABLE 1 ###############

table(dat$convict_type)/sum(!is.na(dat$convict_type), na.rm = T)
ctab <- data.frame(c("Corruption", "Economic Crimes (e.g., fraud) ", "Miscellaneous", "Violence"), as.numeric(table(dat$convict_type)/sum(!is.na(dat$convict_type), na.rm = T)))

print.xtable(xtable(ctab, digits=dig, caption="Types of Charges Brought against Mayors"), type="latex", file="output/table1.tex", include.rownames = F, include.colnames = F, caption.placement = "top")



############# TABLE 2 ###############

md.mayor <- stats::aggregate(cbind(mayor_purge5, regime_independent, appointed) ~ mayor_id, data=dat, mean, na.action=na.omit, subset=inoffice==1)
md.mayor$mayor_purge5 <- ceiling(md.mayor$mayor_purge5)

prop.table(table(ceiling(md.mayor$regime_independent), md.mayor$mayor_purge5), 1)
rtab.num <- cbind(rowSums(table(ceiling(md.mayor$regime_independent), md.mayor$mayor_purge5)), table(ceiling(md.mayor$regime_independent), md.mayor$mayor_purge5)[,2])
rtab <- data.frame(c("Oppositional", "Independent", "Regime-Affiliated"), rtab.num, as.numeric(prop.table(table(ceiling(md.mayor$regime_independent), md.mayor$mayor_purge5), 1)[,2]), stringsAsFactors = F)
colnames(rtab) <- NULL
rownames(rtab) <- NULL

prop.table(table(ceiling(md.mayor$appointed), md.mayor$mayor_purge5), 1)
etab.num <- cbind(rowSums(table(ceiling(md.mayor$appointed), md.mayor$mayor_purge5)), table(ceiling(md.mayor$appointed), md.mayor$mayor_purge5)[,2])
etab <- data.frame(c("Popularly Elected Mayors", "Non-Elected Mayors"), etab.num, as.numeric(prop.table(table(ceiling(md.mayor$appointed), md.mayor$mayor_purge5), 1)[,2]), stringsAsFactors = F)
colnames(etab) <- NULL
rownames(etab) <- NULL

dat.mayor2 <- stats::aggregate(cbind(mayor_purge5, city_manager) ~ mayor_id, data=dat, mean, na.action=na.omit, subset=(inoffice==1 & mayormodel==3))
dat.mayor2$mayor_purge5 <- ceiling(dat.mayor2$mayor_purge5)
dat.mayor2$city_manager <- ceiling(dat.mayor2$city_manager)
dat.mayor3 <- stats::aggregate(cbind(mayor_purge5) ~ city_manager, data=dat.mayor2, mean, na.action=na.omit)
dat.mayor3.num <- stats::aggregate(cbind(mayor_purge5) ~ city_manager, data=dat.mayor2, sum, na.action=na.omit)

ctab.num <- cbind(table(dat.mayor2$city_manager), dat.mayor3.num[,2])
ctab <- data.frame(c("Heads of Cities", "City Managers"), ctab.num, as.numeric(dat.mayor3[,2]), stringsAsFactors = F)
colnames(ctab) <- NULL
rownames(ctab) <- NULL

combotab <- base::rbind(as.matrix(rtab), as.matrix(etab), as.matrix(ctab))
print.xtable(xtable(combotab, digits=dig, caption="Proportion of Mayors Ever Arrested, by Mayor Characteristics"), type="latex", file="output/table2.tex", include.rownames = F, include.colnames = F, caption.placement = "top")


######################################## MAIN RESULTS TABLES #################################################

############# TABLE 3 ###############

t3_m1 <- glm(mayor_purge5 ~ ewi_reg + ewi_force_kgb_prosecutor_mvd + ewi_fed2 + ewi_bus + outsider_dummy + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t3_m1.mfx <- margins::margins(t3_m1, vcov=sandwich::vcovHC(t3_m1, type="HC0", cluster="regionid"))
t3_m1.ses <- summary(t3_m1.mfx)$SE
names(t3_m1.ses) <- str_remove(names(t3_m1.ses), "Var_dydx_")
t3_m1.nummayor <- length(table(dat[rownames(t3_m1$mayormodel),"mayor_id"]))
t3_m1.numcity <- length(table(dat[rownames(t3_m1$mayormodel),"city_id"]))
t3_m1.numregion <- length(table(dat[rownames(t3_m1$mayormodel),"regionid"]))

t3_m2 <- glm(mayor_purge5 ~ rwork1 + rwork2 + rwork3 + rwork4 + rwork5 + outsider_dummy + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t3_m2.mfx <- margins::margins(t3_m2, vcov=sandwich::vcovHC(t3_m2, type="HC0", cluster="regionid"))
t3_m2.ses <- summary(t3_m2.mfx)$SE
names(t3_m2.ses) <- str_remove(names(t3_m2.ses), "Var_dydx_")
t3_m2.nummayor <- length(table(dat[rownames(t3_m2$mayormodel),"mayor_id"]))
t3_m2.numcity <- length(table(dat[rownames(t3_m2$mayormodel),"city_id"]))
t3_m2.numregion <- length(table(dat[rownames(t3_m2$mayormodel),"regionid"]))

t3_m3 <- glm(mayor_purge5 ~ appointed + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t3_m3.mfx <- margins::margins(t3_m3, vcov=sandwich::vcovHC(t3_m3, type="HC0", cluster="regionid"))
t3_m3.ses <- summary(t3_m3.mfx)$SE
names(t3_m3.ses) <- str_remove(names(t3_m3.ses), "Var_dydx_")
t3_m3.nummayor <- length(table(dat[rownames(t3_m3$mayormodel),"mayor_id"]))
t3_m3.numcity <- length(table(dat[rownames(t3_m3$mayormodel),"city_id"]))
t3_m3.numregion <- length(table(dat[rownames(t3_m3$mayormodel),"regionid"]))

t3_m4 <- glm(mayor_purge5 ~ regime_dummy + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t3_m4.mfx <- margins::margins(t3_m4, vcov=sandwich::vcovHC(t3_m4, type="HC0", cluster="regionid"))
t3_m4.ses <- summary(t3_m4.mfx)$SE
names(t3_m4.ses) <- str_remove(names(t3_m4.ses), "Var_dydx_")
t3_m4.nummayor <- length(table(dat[rownames(t3_m4$mayormodel),"mayor_id"]))
t3_m4.numcity <- length(table(dat[rownames(t3_m4$mayormodel),"city_id"]))
t3_m4.numregion <- length(table(dat[rownames(t3_m4$mayormodel),"regionid"]))

t3_m5 <- glm(mayor_purge5 ~ regime_dummy + appointed + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t3_m5.mfx <- margins::margins(t3_m5, vcov=sandwich::vcovHC(t3_m5, type="HC0", cluster="regionid"))
t3_m5.ses <- summary(t3_m5.mfx)$SE
names(t3_m5.ses) <- str_remove(names(t3_m5.ses), "Var_dydx_")
t3_m5.nummayor <- length(table(dat[rownames(t3_m5$mayormodel),"mayor_id"]))
t3_m5.numcity <- length(table(dat[rownames(t3_m5$mayormodel),"city_id"]))
t3_m5.numregion <- length(table(dat[rownames(t3_m5$mayormodel),"regionid"]))


labels_t3 <- c("Ever worked: Regional Gov't", "Ever worked: Siloviki", "Ever worked: Federal Gov't", "Ever worked: Business", "Most recent work: Business", "Most recent work: Local admin", "Most recent work: Regional admin", "Most recent work: Regional legis", "Most recent work: Federal govt", "Outsider", "Non-Elected", "Regime Affiliated", "New Governor", "Percent Ethnic Russian (region)", "Log City Population", "Petrov-Titkov Democracy (region)")

stargazer(list(t3_m1, t3_m2, t3_m3, t3_m4, t3_m5), coef=list(summary(t3_m1.mfx)$AME, summary(t3_m2.mfx)$AME, summary(t3_m3.mfx)$AME, summary(t3_m4.mfx)$AME, summary(t3_m5.mfx)$AME), se=list(t3_m1.ses, t3_m2.ses, t3_m3.ses, t3_m4.ses, t3_m5.ses), type="latex", style="apsr", covariate.labels=labels_t3, font.size="footnotesize", label="arresttable1", dep.var.labels=dv, align=T, no.space=T, order=NULL, add.lines=list(c("Num mayors", t3_m1.nummayor, t3_m2.nummayor, t3_m3.nummayor, t3_m4.nummayor, t3_m5.nummayor), c("Num cities", t3_m1.numcity, t3_m2.numcity, t3_m3.numcity, t3_m4.numcity, t3_m5.numcity), c("Num regions", t3_m1.numregion, t3_m2.numregion, t3_m3.numregion, t3_m4.numregion, t3_m5.numregion)), out="output/table3.tex", notes=note, notes.append=F, omit.stat=omstat, table.placement=tplc, omit=c("factor", "tenure"), column.sep.width = csw, title="Personal Connections, Regime Affiliation, and Mayor Arrest")


############# TABLE 4 ###############

t4_m1 <- glm(mayor_purge5 ~ regime_independent + regionvoteshare_ur + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t4_m1.mfx <- margins::margins(t4_m1, vcov=sandwich::vcovHC(t4_m1, type="HC0", cluster="regionid"))
t4_m1.ses <- summary(t4_m1.mfx)$SE
names(t4_m1.ses) <- str_remove(names(t4_m1.ses), "Var_dydx_")
t4_m1.nummayor <- length(table(dat[rownames(t4_m1$mayormodel),"mayor_id"]))
t4_m1.numcity <- length(table(dat[rownames(t4_m1$mayormodel),"city_id"]))
t4_m1.numregion <- length(table(dat[rownames(t4_m1$mayormodel),"regionid"]))

t4_m2 <- glm(mayor_purge5 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t4_m2.mfx <- margins::margins(t4_m2, vcov=sandwich::vcovHC(t4_m2, type="HC3", cluster="regionid"))
t4_m2.ses <- summary(t4_m2.mfx)$SE
names(t4_m2.ses) <- str_remove(names(t4_m2.ses), "Var_dydx_")
t4_m2.nummayor <- length(table(dat[rownames(t4_m2$mayormodel),"mayor_id"]))
t4_m2.numcity <- length(table(dat[rownames(t4_m2$mayormodel),"city_id"]))
t4_m2.numregion <- length(table(dat[rownames(t4_m2$mayormodel),"regionid"]))

t4_m3 <- glm(mayor_purge5 ~  mayormargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001 & regime_independent==2))
t4_m3.mfx <- margins::margins(t4_m3, vcov=sandwich::vcovHC(t4_m3, type="HC0", cluster="regionid"))
t4_m3.ses <- summary(t4_m3.mfx)$SE
names(t4_m3.ses) <- str_remove(names(t4_m3.ses), "Var_dydx_")
t4_m3.nummayor <- length(table(dat[rownames(t4_m3$mayormodel),"mayor_id"]))
t4_m3.numcity <- length(table(dat[rownames(t4_m3$mayormodel),"city_id"]))
t4_m3.numregion <- length(table(dat[rownames(t4_m3$mayormodel),"regionid"]))

t4_m4 <- glm(mayor_purge5 ~  regime_independent + regionvoteshare_ur + mayormargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
t4_m4.mfx <- margins::margins(t4_m4, vcov=sandwich::vcovHC(t4_m4, type="HC0", cluster="regionid"))
t4_m4.ses <- summary(t4_m4.mfx)$SE
names(t4_m4.ses) <- str_remove(names(t4_m4.ses), "Var_dydx_")
t4_m4.nummayor <- length(table(dat[rownames(t4_m4$mayormodel),"mayor_id"]))
t4_m4.numcity <- length(table(dat[rownames(t4_m4$mayormodel),"city_id"]))
t4_m4.numregion <- length(table(dat[rownames(t4_m4$mayormodel),"regionid"]))

t4_m5 <- glm(mayor_purge5 ~  mayormargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001 & regime_independent<2))
t4_m5.mfx <- margins::margins(t4_m5, vcov=sandwich::vcovHC(t4_m5, type="HC0", cluster="regionid"))
t4_m5.ses <- summary(t4_m5.mfx)$SE
names(t4_m5.ses) <- str_remove(names(t4_m5.ses), "Var_dydx_")
t4_m5.nummayor <- length(table(dat[rownames(t4_m5$mayormodel),"mayor_id"]))
t4_m5.numcity <- length(table(dat[rownames(t4_m5$mayormodel),"city_id"]))
t4_m5.numregion <- length(table(dat[rownames(t4_m5$mayormodel),"regionid"]))


labels_t4 <- c("Regime Affiliated", "UR Vote Share in Regional Elections, by mayor" ,"Mayor Margin of Victory", "New Governor", "Percent Ethnic Russian (region)", "Log City Population", "Petrov-Titkov Democracy (region)")

stargazer(list(t4_m1, t4_m2, t4_m4, t4_m3, t4_m5), coef=list(summary(t4_m1.mfx)$AME, summary(t4_m2.mfx)$AME, summary(t4_m4.mfx)$AME, summary(t4_m3.mfx)$AME, summary(t4_m5.mfx)$AME), se=list(t4_m1.ses, t4_m2.ses, t4_m4.ses, t4_m3.ses, t4_m5.ses), type="latex", style="apsr", covariate.labels=labels_t4, font.size="footnotesize", label="arresttable2", dep.var.labels=dv, align=T, no.space=T, order=NULL, add.lines=list(c("Num mayors", t4_m1.nummayor, t4_m2.nummayor, t4_m4.nummayor, t4_m3.nummayor, t4_m5.nummayor), c("Num cities", t4_m1.numcity, t4_m2.numcity, t4_m4.numcity, t4_m3.numcity, t4_m5.numcity), c("Num regions", t4_m1.numregion, t4_m2.numregion, t4_m4.numregion, t4_m3.numregion, t4_m5.numregion)), out="output/table4.tex", notes=note, notes.append=F, omit.stat=omstat, table.placement=tplc, omit=c("factor", "tenure"), column.sep.width = csw, title="Regime Service, Mayor Popularity, and Mayor Arrest")


############# TABLE 5 ###############


t5_m1 <- glm(mayor_purge5 ~  mayormargin + factor(year) + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t5_m1.mfx <- margins::margins(t5_m1, vcov=sandwich::vcovHC(t5_m1, type="HC0", cluster="regionid"))
t5_m1.ses <- summary(t5_m1.mfx)$SE
names(t5_m1.ses) <- str_remove(names(t5_m1.ses), "Var_dydx_")
t5_m1.nummayor <- length(table(dat[rownames(t5_m1$mayormodel),"mayor_id"]))
t5_m1.numcity <- length(table(dat[rownames(t5_m1$mayormodel),"city_id"]))
t5_m1.numregion <- length(table(dat[rownames(t5_m1$mayormodel),"regionid"]))

t5_m2 <- glm(mayor_purge5 ~ regime_independent + mayormargin + factor(year) + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t5_m2.mfx <- margins::margins(t5_m2, vcov=sandwich::vcovHC(t5_m2, type="HC0", cluster="regionid"))
t5_m2.ses <- summary(t5_m2.mfx)$SE
names(t5_m2.ses) <- str_remove(names(t5_m2.ses), "Var_dydx_")
t5_m2.nummayor <- length(table(dat[rownames(t5_m2$mayormodel),"mayor_id"]))
t5_m2.numcity <- length(table(dat[rownames(t5_m2$mayormodel),"city_id"]))
t5_m2.numregion <- length(table(dat[rownames(t5_m2$mayormodel),"regionid"]))


t5_m3 <- glm(mayor_purge5 ~ regime_independent + pop_margin + mayormargin + logpop + factor(year) + newgovernor + pctrussian + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
t5_m3.mfx <- margins::margins(t5_m3, vcov=sandwich::vcovHC(t5_m3, type="HC0", cluster="regionid"))
t5_m3.ses <- summary(t5_m3.mfx)$SE
names(t5_m3.ses) <- str_remove(names(t5_m3.ses), "Var_dydx_")
t5_m3.nummayor <- length(table(dat[rownames(t5_m3$mayormodel),"mayor_id"]))
t5_m3.numcity <- length(table(dat[rownames(t5_m3$mayormodel),"city_id"]))
t5_m3.numregion <- length(table(dat[rownames(t5_m3$mayormodel),"regionid"]))


t5_m4 <- glm(mayor_purge5 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001 & regionid!=10))
t5_m4.mfx <- margins::margins(t5_m4, vcov=sandwich::vcovHC(t5_m4, type="HC0", cluster="regionid"))
t5_m4.ses <- summary(t5_m4.mfx)$SE
names(t5_m4.ses) <- str_remove(names(t5_m4.ses), "Var_dydx_")
t5_m4.nummayor <- length(table(dat[rownames(t5_m4$mayormodel),"mayor_id"]))
t5_m4.numcity <- length(table(dat[rownames(t5_m4$mayormodel),"city_id"]))
t5_m4.numregion <- length(table(dat[rownames(t5_m4$mayormodel),"regionid"]))


t5_m5 <- glm(mayor_purge ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
t5_m5.mfx <- margins::margins(t5_m5, vcov=sandwich::vcovHC(t5_m5, type="HC0", cluster="regionid"))
t5_m5.ses <- summary(t5_m5.mfx)$SE
names(t5_m5.ses) <- str_remove(names(t5_m5.ses), "Var_dydx_")
t5_m5.nummayor <- length(table(dat[rownames(t5_m5$mayormodel),"mayor_id"]))
t5_m5.numcity <- length(table(dat[rownames(t5_m5$mayormodel),"city_id"]))
t5_m5.numregion <- length(table(dat[rownames(t5_m5$mayormodel),"regionid"]))

t5_m6 <- glm(mayor_purge6 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
t5_m6.mfx <- margins::margins(t5_m6, vcov=sandwich::vcovHC(t5_m6, type="HC0", cluster="regionid"))
t5_m6.ses <- summary(t5_m6.mfx)$SE
names(t5_m6.ses) <- str_remove(names(t5_m6.ses), "Var_dydx_")
t5_m6.nummayor <- length(table(dat[rownames(t5_m6$mayormodel),"mayor_id"]))
t5_m6.numcity <- length(table(dat[rownames(t5_m6$mayormodel),"city_id"]))
t5_m6.numregion <- length(table(dat[rownames(t5_m6$mayormodel),"regionid"]))



labels_t5 <- c("Regime Affiliated", "City Population x Mayor Margin of Victory", "Mayor Margin of Victory", "Percent Ethnic Russian (region)", "Log City Population", "Petrov-Titkov Democracy (region)")

stargazer(list(t5_m1, t5_m2, t5_m4, t5_m3, t5_m6, t5_m5), coef=list(summary(t5_m1.mfx)$AME, summary(t5_m2.mfx)$AME, summary(t5_m4.mfx)$AME, summary(t5_m3.mfx)$AME, summary(t5_m6.mfx)$AME, summary(t5_m5.mfx)$AME), se=list(t5_m1.ses, t5_m2.ses, t5_m4.ses, t5_m3.ses, t5_m6.ses, t5_m5.ses), type="latex", style="apsr", covariate.labels=labels_t5, font.size="footnotesize", label="robust1", dep.var.labels=dv, align=T, no.space=T, order=NULL, add.lines=list(c("Num mayors", t5_m1.nummayor, t5_m2.nummayor, t5_m4.nummayor, t5_m3.nummayor, t5_m6.nummayor, t5_m5.nummayor), c("Num cities", t5_m1.numcity, t5_m2.numcity, t5_m4.numcity, t5_m3.numcity, t5_m6.numcity, t5_m5.numcity), c("Num regions", t5_m1.numregion, t5_m2.numregion, t5_m4.numregion, t5_m3.numregion, t5_m6.numregion, t5_m5.numregion)), out="output/table5.tex", notes=note, notes.append=F, omit.stat=omstat, table.placement=tplc, omit=c("year", "tenure"), column.sep.width = csw, title="Robustness")


######################################## APPENDIX TABLES #################################################

############# TABLE A1 ###############


ta1_m1 <- glm(mayor_purge1 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
ta1_m1.mfx <- margins::margins(ta1_m1, vcov=sandwich::vcovHC(ta1_m1, type="HC0", cluster="regionid"))
ta1_m1.ses <- summary(ta1_m1.mfx)$SE
names(ta1_m1.ses) <- str_remove(names(ta1_m1.ses), "Var_dydx_")
ta1_m1.nummayor <- length(table(dat[rownames(ta1_m1$mayormodel),"mayor_id"]))
ta1_m1.numcity <- length(table(dat[rownames(ta1_m1$mayormodel),"city_id"]))
ta1_m1.numregion <- length(table(dat[rownames(ta1_m1$mayormodel),"regionid"]))

ta1_m2 <- glm(mayor_purge3 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
ta1_m2.mfx <- margins::margins(ta1_m2, vcov=sandwich::vcovHC(ta1_m2, type="HC0", cluster="regionid"))
ta1_m2.ses <- summary(ta1_m2.mfx)$SE
names(ta1_m2.ses) <- str_remove(names(ta1_m2.ses), "Var_dydx_")
ta1_m2.nummayor <- length(table(dat[rownames(ta1_m2$mayormodel),"mayor_id"]))
ta1_m2.numcity <- length(table(dat[rownames(ta1_m2$mayormodel),"city_id"]))
ta1_m2.numregion <- length(table(dat[rownames(ta1_m2$mayormodel),"regionid"]))

ta1_m3 <- glm(mayor_purge4 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
ta1_m3.mfx <- margins::margins(ta1_m3, vcov=sandwich::vcovHC(ta1_m3, type="HC0", cluster="regionid"))
ta1_m3.ses <- summary(ta1_m3.mfx)$SE
names(ta1_m3.ses) <- str_remove(names(ta1_m3.ses), "Var_dydx_")
ta1_m3.nummayor <- length(table(dat[rownames(ta1_m3$mayormodel),"mayor_id"]))
ta1_m3.numcity <- length(table(dat[rownames(ta1_m3$mayormodel),"city_id"]))
ta1_m3.numregion <- length(table(dat[rownames(ta1_m3$mayormodel),"regionid"]))

ta1_m4 <- glm(mayor_purge5 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + dem + logpop + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
ta1_m4.mfx <- margins::margins(ta1_m4, vcov=sandwich::vcovHC(ta1_m4, type="HC0", cluster="regionid"))
ta1_m4.ses <- summary(ta1_m4.mfx)$SE
names(ta1_m4.ses) <- str_remove(names(ta1_m4.ses), "Var_dydx_")
ta1_m4.nummayor <- length(table(dat[rownames(ta1_m4$mayormodel),"mayor_id"]))
ta1_m4.numcity <- length(table(dat[rownames(ta1_m4$mayormodel),"city_id"]))
ta1_m4.numregion <- length(table(dat[rownames(ta1_m4$mayormodel),"regionid"]))

labels_a1 <- c("Regime Affiliated", "Mayor Margin of Victory", "New Governor", "Percent Ethnic Russian (region)", "Log City Population", "Petrov-Titkov Democracy (region)")

stargazer(list(ta1_m4, ta1_m3, ta1_m2, ta1_m1), coef=list(summary(ta1_m4.mfx)$AME, summary(ta1_m3.mfx)$AME, summary(ta1_m2.mfx)$AME, summary(ta1_m1.mfx)$AME), se=list(ta1_m4.ses, ta1_m3.ses, ta1_m2.ses, ta1_m1.ses), type="latex", style="apsr", covariate.labels=labels_a1, font.size="footnotesize", label="robust2_dv", dep.var.labels=dv, align=T, no.space=T, order=NULL, add.lines=list(c("Num mayors", ta1_m4.nummayor, ta1_m3.nummayor, ta1_m2.nummayor, ta1_m1.nummayor), c("Num cities", ta1_m4.numcity, ta1_m3.numcity, ta1_m2.numcity, ta1_m1.numcity), c("Num regions", ta1_m4.numregion, ta1_m3.numregion, ta1_m2.numregion, ta1_m1.numregion)), out="output/table_a1.tex", notes=note, notes.append=F, omit.stat=omstat, table.placement=tplc, omit=c("factor", "tenure"), column.sep.width = csw, title="Additional Robustness Checks")



############# TABLE A2 ###############

ta2_m1 <- glm(mayor_purge5 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001 & year<2011))
ta2_m1.mfx <- margins::margins(ta2_m1, vcov=sandwich::vcovHC(ta2_m1, type="HC0", cluster="regionid"))
ta2_m1.ses <- summary(ta2_m1.mfx)$SE
names(ta2_m1.ses) <- str_remove(names(ta2_m1.ses), "Var_dydx_")
ta2_m1.nummayor <- length(table(dat[rownames(ta2_m1$mayormodel),"mayor_id"]))
ta2_m1.numcity <- length(table(dat[rownames(ta2_m1$mayormodel),"city_id"]))
ta2_m1.numregion <- length(table(dat[rownames(ta2_m1$mayormodel),"regionid"]))

ta2_m2 <- glm(mayor_purge5 ~ regime_independent + mayormargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2010))
ta2_m2.mfx <- margins::margins(ta2_m2, vcov=sandwich::vcovHC(ta2_m2, type="HC0", cluster="regionid"))
ta2_m2.ses <- summary(ta2_m2.mfx)$SE
names(ta2_m2.ses) <- str_remove(names(ta2_m2.ses), "Var_dydx_")
ta2_m2.nummayor <- length(table(dat[rownames(ta2_m2$mayormodel),"mayor_id"]))
ta2_m2.numcity <- length(table(dat[rownames(ta2_m2$mayormodel),"city_id"]))
ta2_m2.numregion <- length(table(dat[rownames(ta2_m2$mayormodel),"regionid"]))

ta2_m3 <- glm(mayor_purge5 ~ regime_independent + mayormargin + regionvoteshare_ur + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu + reg_grow, data=dat, family=binomial, subset=year>2001)
ta2_m3.mfx <- margins::margins(ta2_m3, vcov=sandwich::vcovHC(ta2_m3, type="HC0", cluster="regionid"))
ta2_m3.ses <- summary(ta2_m3.mfx)$SE
names(ta2_m3.ses) <- str_remove(names(ta2_m3.ses), "Var_dydx_")
ta2_m3.nummayor <- length(table(dat[rownames(ta2_m3$mayormodel),"mayor_id"]))
ta2_m3.numcity <- length(table(dat[rownames(ta2_m3$mayormodel),"city_id"]))
ta2_m3.numregion <- length(table(dat[rownames(ta2_m3$mayormodel),"regionid"]))

ta2_m4 <- glm(mayor_purge5 ~ regime_independent + mayormargin + regionvoteshare_ur + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu + city_unemp, data=dat, family=binomial, subset=year>2001)
ta2_m4.mfx <- margins::margins(ta2_m4, vcov=sandwich::vcovHC(ta2_m4, type="HC0", cluster="regionid"))
ta2_m4.ses <- summary(ta2_m4.mfx)$SE
names(ta2_m4.ses) <- str_remove(names(ta2_m4.ses), "Var_dydx_")
ta2_m4.nummayor <- length(table(dat[rownames(ta2_m4$mayormodel),"mayor_id"]))
ta2_m4.numcity <- length(table(dat[rownames(ta2_m4$mayormodel),"city_id"]))
ta2_m4.numregion <- length(table(dat[rownames(ta2_m4$mayormodel),"regionid"]))

ta2_m5 <- glm(mayor_purge5 ~ regime_independent + mayormargin + ds_perc_rejected + regionvoteshare_ur + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
ta2_m5.mfx <- margins::margins(ta2_m5, vcov=sandwich::vcovHC(ta2_m5, type="HC0", cluster="regionid"))
ta2_m5.ses <- summary(ta2_m5.mfx)$SE
names(ta2_m5.ses) <- str_remove(names(ta2_m5.ses), "Var_dydx_")
ta2_m5.nummayor <- length(table(dat[rownames(ta2_m5$mayormodel),"mayor_id"]))
ta2_m5.numcity <- length(table(dat[rownames(ta2_m5$mayormodel),"city_id"]))
ta2_m5.numregion <- length(table(dat[rownames(ta2_m5$mayormodel),"regionid"]))

ta2_m6 <- glm(mayor_purge5 ~ regime_independent + mayormargin + cc_urpct + regionvoteshare_ur + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=year>2001)
ta2_m6.mfx <- margins::margins(ta2_m6, vcov=sandwich::vcovHC(ta2_m6, type="HC0", cluster="regionid"))
ta2_m6.ses <- summary(ta2_m6.mfx)$SE
names(ta2_m6.ses) <- str_remove(names(ta2_m6.ses), "Var_dydx_")
ta2_m6.nummayor <- length(table(dat[rownames(ta2_m6$mayormodel),"mayor_id"]))
ta2_m6.numcity <- length(table(dat[rownames(ta2_m6$mayormodel),"city_id"]))
ta2_m6.numregion <- length(table(dat[rownames(ta2_m6$mayormodel),"regionid"]))


labels_a2 <- c("Regime Affiliated","Mayor Margin of Victory", "UR Pct City Council Vote", "Pct Candidates Rejected", "UR Vote Share in Regional Elections, by mayor", "New Governor", "Percent Ethnic Russian (region)", "Log City Population", "Petrov-Titkov Democracy (region)", "Regional Econ Growth", "City Unemployment Level")

stargazer(list(ta2_m1, ta2_m2, ta2_m3, ta2_m4, ta2_m6, ta2_m5), coef=list(summary(ta2_m1.mfx)$AME, summary(ta2_m2.mfx)$AME, summary(ta2_m3.mfx)$AME, summary(ta2_m4.mfx)$AME, summary(ta2_m6.mfx)$AME, summary(ta2_m5.mfx)$AME), se=list(ta2_m1.ses, ta2_m2.ses, ta2_m3.ses, ta2_m4.ses, ta2_m6.ses, ta2_m5.ses), type="latex", style="apsr", covariate.labels=labels_a2, font.size="footnotesize", label="rr1", dep.var.labels=dv, align=T, no.space=T, order=NULL, add.lines=list(c("Num mayors", ta2_m1.nummayor, ta2_m2.nummayor, ta2_m3.nummayor, ta2_m4.nummayor, ta2_m6.nummayor, ta2_m5.nummayor), c("Num cities", ta2_m1.numcity, ta2_m2.numcity, ta2_m3.numcity, ta2_m4.numcity, ta2_m6.numcity, ta2_m5.numcity), c("Num regions", ta2_m1.numregion, ta2_m2.numregion, ta2_m3.numregion, ta2_m4.numregion, ta2_m6.numregion, ta2_m5.numregion)), out="output/table_a2.tex", notes=note, notes.append=F, omit.stat=omstat, table.placement=tplc, omit=c("factor", "tenure"), column.sep.width = csw, title="Robustness: Effects over Time, Economic Performance, City Council")



############# TABLE A3 ###############

ta3_m1 <- glm(mayor_purge5 ~  putin_popularity + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
ta3_m1.mfx <- margins::margins(ta3_m1, vcov=sandwich::vcovHC(ta3_m1, type="HC0", cluster="regionid"))
ta3_m1.ses <- summary(ta3_m1.mfx)$SE
names(ta3_m1.ses) <- str_remove(names(ta3_m1.ses), "Var_dydx_")
ta3_m1.nummayor <- length(table(dat[rownames(ta3_m1$mayormodel),"mayor_id"]))
ta3_m1.numcity <- length(table(dat[rownames(ta3_m1$mayormodel),"city_id"]))
ta3_m1.numregion <- length(table(dat[rownames(ta3_m1$mayormodel),"regionid"]))


ta3_m2 <- glm(mayor_purge5 ~  mayormargin + putin_popularity + putinXmargin + factor(year) + newgovernor + pctrussian + logpop + dem + mayortenure + mayortenure_sq + mayortenure_cu, data=dat, family=binomial, subset=(year>2001))
ta3_m2.mfx <- margins::margins(ta3_m2, vcov=sandwich::vcovHC(ta3_m2, type="HC0", cluster="regionid"))
ta3_m2.ses <- summary(ta3_m2.mfx)$SE
names(ta3_m2.ses) <- str_remove(names(ta3_m2.ses), "Var_dydx_")
ta3_m2.nummayor <- length(table(dat[rownames(ta3_m2$mayormodel),"mayor_id"]))
ta3_m2.numcity <- length(table(dat[rownames(ta3_m2$mayormodel),"city_id"]))
ta3_m2.numregion <- length(table(dat[rownames(ta3_m2$mayormodel),"regionid"]))


labels_a3 <- c("Mayor Margin of Victory", "Putin Natl Popularity", "Putin Pop x Mayor Margin", "New Governor", "Percent Ethnic Russian (region)", "Log City Population", "Petrov-Titkov Democracy (region)")

stargazer(list(ta3_m1, ta3_m2), coef=list(summary(ta3_m1.mfx)$AME, summary(ta3_m2.mfx)$AME), se=list(ta3_m1.ses, ta3_m2.ses), type="latex", style="apsr", covariate.labels=labels_a3, font.size="footnotesize", label="putinpop", dep.var.labels=dv, align=T, no.space=T, order=NULL, add.lines=list(c("Num mayors", ta3_m1.nummayor, ta3_m2.nummayor), c("Num cities", ta3_m1.numcity, ta3_m2.numcity), c("Num regions", ta3_m1.numregion, ta3_m2.numregion)), out="output/table_a3.tex", notes=note, notes.append=F, omit.stat=omstat, table.placement=tplc, omit=c("factor", "tenure"), column.sep.width = csw, title="Putin Popularity and Mayor Arrest")


